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Abstract 

Healthy and sick individuals {A and B particles) diffuse independently 
with diffusion constants Da and Db- Sick individuals upon encounter 
infect healthy ones (at rate fc), but may also spontaneously recover (at 
rate 1/t). The propagation of the epidemic therefore couples to the 
fluctuations in the total population density. Global extinction occurs 
below a critical value pc of the spatially averaged total density. The 
epidemic evolves as the diffusion-reaction-decay process A + B ~> 
2B, B ^ A, for which we write down the field theory. The stationary 
state properties of this theory when Da = Db were obtained by Kree 
et al. The critical behavior for Da < Db is governed by a new fixed 
point. We calculate the critical exponents of the stationary state in an 
e expansion, carried out by Wilson renormalization, below the critical 
dimension dc = 4. We then go on to to obtain the critical initial 
time behavior at the extinction threshold, both for Da = Db and 
Da < Db- There is nonuniversal dependence on the initial particle 
distribution. The case Da > Db remains unsolved. 
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1 Introduction and phenomenological analysis 



1.1 Introduction 

Wc study the propagation of an epidemic in a population of fluctuating den- 
sity, with parameters set such that the epidemic is close to becoming extinct. 
The model is as follows. Healthy and sick individuals, also called A and B 
particles, diffuse freely and independently on a lattice of dimension d. Their 
diffusion constants Da and Db are in general unequal. Sick individuals upon 
encounter infect healthy ones (at a rate k), but may also spontaneously re- 
cover (at a rate 1/r) . The epidemic therefore tends to propagate in regions 
of high total population density and to become locally extinct in regions 
of low density. Below a critical density pc global extinction occurs, and an 
important question is how the epidemic evolves near this threshold, which 
is analogous to a phase transition point in equilibrium critical phenomena. 
The role of the order parameter is played by the B particle density pB- 
In particular, at threshold strong local fluctuations arc expected, with a 
potentially decisive impact on the global critical behavior. 

The evolution of the epidemic may be formulated as the competition 
between a reaction and a decay process of two diffusing chemical species, 

A + B^2B, B^-^A (1.1) 

in which A and B stand for a healthy and a sick individual, respectively. 

Let rrii and rij be the number of A and B particles, respectively, that 
occupy site i. We write m = {nii} and n = {rii} for the two sets of oc- 
cupation numbers and denote by P{m, n; t) the probability that at time t 
the system is characterized by exactly these sets. The rules of the process 
formulated above then define a master equation for the time evolution of 
P{m, n; t) with fully specified transition rates. The problem depends on the 
parameters Da, Db, k, 1/r, and p. After scaling one may set A; = 1. 

This problem is but one instance of a large class of such processes that 
all have the following common features. They are initially formulated in 
terms of a master equation which is hard, if not impossible, to solve. The 
traditional approach is to write down rate equations for the space and time 
dependent densities of the various species. This leads to numerous interest- 
ing systems of partial differential equations (reaction-diffusion equations) 
that continue to be subjects of investigation. 

The assumptions that lead to such reaction-diffusion equations can be 
justified only in sufficiently high spatial dimension d. The reason is that 
these equations constitute a mean field description: the reaction terms in- 
volve the products of the average densities, rather than the average of the 
product of the fluctuating densities. 

The general problem of statistical mechanics is therefore to determine 
the critical dimension d = d^ below which the mean-field approach becomes 
qualitatively incorrect, and to describe the modifications that then appear. 
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Recent progress on such problems is due to the appUcation of the renormal- 
ization group (RG) to the formal solution of the master equation, rewritten 
as a path integral and cast in the form of a field theory. This will also be our 
procedure. It allows one, in an expansion in the small parameter dc — d = e, 
to study the scale invariant regimes of the problem and determine the asso- 
ciated exponents. 

For the particular process of Eq. (|L1D the upper critical dimension ap- 
pears to be dc = 4. Questions of interest are first of all about the properties 
of the stationary state and secondly about how the system, prepared in a 
specified initial state, relaxes towards stationarity. Of particular interest 
is the initial state with the A and B particles distributed randomly with 
densities and p^. In that case, and in the limit of an infinite lattice, the 
nii and nj are independent Poissonian random variables with averages p^ 
and p^, respectively. 

We now relate our work to existing literature. We shall arrive at a field 
theory described by an action S (see Eq. ( 1.16D ) in which 



is a key parameter. The stationary state of a special case - viz. with 
Da = Db ~ was studied recently by Kree, Schaub, and Schmittmann [||]. 
These authors' physical interpretation is different: They start from a pair 
of Langevin equations postulated for the evolution of a population density 
in a polluted environment. The interest of both this work and of Ref. j^] 
lies in the coupling between the fluctuations of the total density and those 
of the density of the species that is trying to survive (in our case the B 
particles). If the fluctuations in the total density are suppressed, the full 
action S becomes independent of p. and reduces to an action 5dp known to 
represent both the problem of directed percolation and the Schlogl reaction. 
^DP has flrst appeared in particle physics as Reggeon fleld theory. 

This work goes beyond the existing results in the following ways. 

(i) The ratio p appears to be relevant in the sense of the RG. For p < 
the critical behavior is governed by a new flxed point. We determine the 
new exponents to leading order in e = 4 — d. For p > the RG equations 
do not have a flxed point and we speculate about the interpretation of this 
fact. 

(a) We consider relaxation to the stationary state at the critical density 
p = Pc first for the action ^dp of directed percolation and then for the full 
action S both with p = and p < 0. For asymptotically long times the order 
parameter psit) relaxes towards zero as a power law of time whose expo- 
nent is a simple combination {f3/vz in the usual notation) of the stationary 
state exponents. There is, however, an initial time regime - which in the 
limit /Og i becomes arbitrarily long - where the order parameter increases 
as , where 0' is a new and independent critical initial slip exponent. We 
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calculate this exponent for an initial state in which the and rij are inde- 
pendent but not necessarily Poissonian random variables. We show that 0' 
is nonuniversal, as it depends continuously on the width of the distribution. 
The concept of critical initial slip was introduced by Janssen, Schaub, and 
Schmittmann |Q] (see also Janssen Q) and it is here for the first time, to 
our knowledge, that such an exponent is calculated for a reaction-diffusion 
system. 

(in) Finally, we innovate with respect to the existing literature by not 
following the usual field-theoretic methods. Instead, we perform a Wilson 
renormalization, that is, we successively integrate out shells of large momen- 
tum. We employ the dynamical version of Wilson's method developed by 
Hohenberg and Halperin [^. This method does not require any advanced 
knowledge of field theory. In the course of our work we check that in known 
special cases earlier results are reproduced. 



1.2 Mean-field theory 

The mean field equations for the densities pA and pB are 

— — = Da^Pa + -PB - k pAPB 
ot T 

= Db^PB - -PB + k PAPB (1-3) 

By summing them one sees that the total density p = pA + Pb conserved. 
The spatially homogeneous solutions of these equations are easily found. 
The total density p appears to have a threshold value p = p™^ = (kr)^^ 
below which the stationary state consists of only A particles. Above this 
threshold density the all-^ state is still stationary but unstable, and the 
stable stationary state has a nonzero density of B particles. Explicitly, 
Eq. (O) leads to 



PAioo)=pf^, pB{oo)=p-pf i0Tp>pf 

Pa{oo) = p, /)b(oo)=0 ior p < p'f 



(1.4) 



If we write, in analogy with critical phenomena in thermodynamic systems, 

PB{^)''^''ip-pcf (1.5) 

then one sees that mean field theory gives (3=1. The relaxation to the sta- 
tionary state at criticality from a homogeneous initial state with i?-particle 
density is easily derived from Eq. (1.3), 



PB{t) = Y^hnz (1-6) 



Again in analogy with critical relaxation in thermodynamic systems we de- 
fine the exponents zv and 9' by 

PB{t) ~ t-PI^^^^ for t » {p%k)-^ (1.7) 
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PB{t) ~ Pst^' for p%kt < 1 (1.8) 

Here 9' is the critical initial slip exponent introduced by Janssen et al. [||, ^. 
Eq. ( |1.6|) shows that the mean field values are zu = 1 and 9' = 0. In 
the following sections, after showing that the reaction-diffusion process of 
Eq. ( |1 . iD has upper critical dimension dc = 4, we will compute corrections to 
the mean field values of these and other exponents to first order in e = 4 — d. 



1.3 The master equation solved in terms of a path integral 

In order to go beyond mean field theory, one has to start from the master 
equation defined implicitly in subsection ^T^. It is convenient to introduce an 
orthonormal basis of states \m,n) and work in the Hilbert space generated 
by this basis. The probability distribution P{m, n; t) is then associated with 
the state 

\P{t)) = J2 P{rn, n- t)\m, n) (1.9) 



For an appropriately constructed time evolution operator H the master 
equation for P{m, n; t) is equivalent to the imaginary time Schrodinger equa- 
tion 

^|P(t)) = -H\P{t)) (1.10) 
From the formal solution of Eq. ( |1.10| ), 

|P(t)) = e-^*|P(0)) (1.11) 



all properties of the reaction-diffusion system may be derived. Eq. ( 1.11 ) de- 
termines in particular in the limit t — > oo the stationary state. The combined 
efforts of many workers, e.g. ^, 0], have led to methods of converting the 



formal solution of type ( 1.11] ) into expressions that can be analyzed. We 



shall use here one such method, which has now become common. The first 
step is to express H in terms of creation and annihilation operators for the 
two particle species: a| and Oj for the A particles, and bj and 6j for the B 
particles. These operators can be defined such that they have boson commu- 



tation relations. Next, the exponential in ( 1.11 ) is subjected to time slicing 
with the aid of a coherent state representation of the harmonic oscillator 
operators. As a result, and in the limit where the site index becomes a 
continuous space vector x, there appear space and time dependent classical 
fields a*(x,t), a{x,t) and b*{x,t), b{x,t) that are directly associated with 



the two particle types. The solution (|1.11 ) finally takes the form of a path 



integral on these fields, weighted with the exponential of an action S: 

\P{t)) = J VaVa*VbVb* Q-S[a,a* ,br\\p^Q^^^^ ^^^^2) 

The exact way in which such physical quantities as local densities, correlation 
and response functions can be obtained from the fields is detailed in Q. 
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1.4 Action and scaling dimensionalities 

The full action for the problem of Eq. is 



S[a,a*,b, b*] 



d'^xdt 



a*{dt - DAA)a + b*{dt - £'bA)6 



+ kab{a* - b*)b* - -{a* -b*)b 



pla*{x, 0) + pW{x, 0) + a{x, T) + b{x, T) 



1.13) 



Here p\ and are the initial A and B particle densities, and Poissonian 
distribution of the initial occupation numbers has been assumed; and the 
time integral runs through a fixed interval [0,T]. The averages with weight 
e""^ of the fields a{x,t) and b{x,t) are equal fS], pip] to the local densities 
of the A and B particles. In order to cast Eq. (|1.13 ) in a form suitable to 



subsequent analysis, we first change from a* and b* to the variables 

a = a* -1 and b = b* - 1 (1.14) 

the effect of which is to get rid of the boundary terms at time t = T. Then 
we introduce the fields 



a + b 
b 



pQ{t) 



a 
b- 



(1.15) 



where the parameter p denotes as before the average total particle density, 
and Q{t) is the unit step function. After further rescaling of the fields and 
time one finds the expression for the action that will be the starting point 
in the following sections. 



s[ip,ip,i',tp] 



d'^xdt Tp{dt - A)ip + i;{dt + A(cj - A))V' 

+ /i^A-i/; + gipipiip — V') + uiljip{ip + Tp) 

+ f 1(^-0)^ + V2ipip{ipTp - ipip) + vsipTpiptp (1.16) 



Vl = 



g = k^/DA 

\a = kipf - p)/D, 



u 
A 



-k^p/DA 
Db/Da 



:i.i7) 



In terms of the original parameters of the master equation the coupling 
constants of this action are given by 

1 - Db/Da 

-- V2 = -v-s = k/DA 
Pb = P%/^/P 

If one omits from the action Eq. ( |1.16| ) the initial time term proportional 
to then the remainder is, for = 0, invariant under the time reversal 
symmetry 



(1.18) 



;0(x,t) - 


^ -tlj{x, 


-t) 


'4){x,t) - 


^ -tlj{x, 


-t) 


ip{x,t) - 


^ ^(x. 


-t) 


W{x,t) - 


^ (/?(x. 


-t) 
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The breaking of this symmetry for 7^ 0, that is, when the diffusion con- 
stants Da and Db are different, wih be seen to have far reaching conse- 
quences for the critical behavior of this system. The naive scahng dimen- 
sionahties of the coupHng constants are, in powers of an inverse length, 

[g] = [u] = 2-d/2, N = 2-d, [H=0 (1.19) 

This shows that dc = 4 is the upper critical dimension for this problem and 
that in an expansion around d = 4 we may drop the Vi vertices. 

The terms occurring in Eq. (|1.16| ) may be given the following physical inter- 
pretation. The g vertex, which involves only the B fields %[) and tp, accounts 
for fluctuations in the contamination process A + B ^ 2B due to fluctu- 
ations in the B particle density, and as though the A particles formed a 
homogeneous background. The sum ^p + 'ip may be shown to be equal to the 
fluctuating part of the total density. The u vertex couples this to the density 
of B particles, expressing that the reaction A + B ^ 2B is also affected by 
fluctuations of the A particle density. 



1.5 A reduced problem: Directed Percolation 

If in the action S of Eq. ( |1.16| ) we replace the fields if and Tp by their mean 
field values if = and ^ = 0, and drop the fourth order terms that were 
argued above to be irrelevant, the result is the action 5dp known in particle 
physics as the Reggeon field theory and given by 



(1 



.20) 



Physically this means that we have constrained the total density to be con- 
stant (equal to p) in space and time. If the B density is much lower than the 
A density, as is the case near criticality, this constraint means that the B 
particles interact with a continuous background. The first link between the 
action S'dp and statistical mechanics was established by Cardy and Sugar 
[Q, who showed that it describes the problem of directed bond percola- 
tion. Subsequently it was shown by Grassberger and Sundermeyer Q and 
by Grassberger and de la Torre [10| (see also Janssen |jll|]) that this action 



also describes Schlogl's chemical reaction [12]. The e expansion for the ex- 



ponents of the action S'dp was known from field theory (||T3] and references 



therein) and has been extended by Janssen |11]. 

The remainder of this article is organized as follows. In Section ^, for a 
twofold reason, we return to the action S'dp- We consider for the first time its 
relaxation from an initial state and calculate the initial critical slip exponent 
6' . Secondly, since throughout this work we employ time dependent Wilson 
renormalization, we explain our method on the example of the relatively 
simple action Sdp, after which we apply it with only few extra comments 



to the full action S of Eq. (|1.16D 
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In section y we return to the full action S. We determine its stationary 
state properties, which are then used to study the relaxation from the initial 
state. We also study the effect of non-Poissonian initial particle distribu- 
tions. 

2 Directed Percolation: Wilson renormalization 
and critical initial slip exponent 

The stationary state of the action Eq. ( |1.20| ) exhibits critical behavior that 
has been studied by several authors both within the framework of an e- 
expansion around its upper critical dimension d = 4. The analysis of the ac- 
tion Eq. ( 1.20 ) is usually performed {e.g. fs], |ll|) by means of field-theoretic 



methods, based upon the perturbative expansion of vertex functions and 
dimensional regularization, and leading to Callan-Symanzik equations. In 
this work, instead, we apply Wilson's dynamic renormalization group, that 
is, successively eliminate short wavelength degrees of freedom. This method 
has already been used successfully to study the dynamics of spin systems for 
instance; it is applied here for the first time to a reaction-diffusion problem. 

2.1 Wilson's dynamic renormalization group 

The extension of Wilson's renormalization procedure to time dependent sys- 
tems has been described by Hohenberg and Halperin Q. We recall here how 



this method works on the example of the action ^dp of Eq. (1.20). For a 
study restricted to the stationary state the boundary term at t = may be 
suppressed. The action may then be written as 5dp = 5o -|- Smt with the 
free term Sq given by 



,So[V,V] = J Hdt + X{a-A))^ (2.1) 
The free propagator, in terms of the spatial Fourier components, is 

ii^i-k, t')i;{k, t)) = e{t - f')e-('='+-)(*-*') (2.2) 

One is interested in averages carried out with the weight e"'^^^ . In Wilson's 
procedure this weight is first integrated over the Fourier components of the 
fields with wavevectors A: in a momentum shell 0,\ defined by 

nA = {k\A/b < k < A} (2.3) 

We denote these fields for brevity by ip^ and '0> and the remaining fields by 
'(/'< and '(/'<• The scale factor b will be taken infinitesimally close to 1, and 
the momentum space cut-off A is of order unity. One obtains an effective 
action Sf^ 

g-S£ff [,^<,?<] ^ r p[^^,^^]e-^DP (2.4) 
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In practice the integration step is carried out perturbatively in Si^t by using 
the cumulant expansion 

exp i-S^l) = exp (-(5int)' + ^{Sfj', - ^{Sfj', + ...) 

^- (2.5) 

The index c indicates a cumulant and (...)' stands for an average on the > 
fields with weight e~^°. Therefore there remain only Gaussian integrations 
to perform. 

Wick's theorem ensures that the integrations reduce to pair contractions, 
and the graphical rules for dealing with these are easily found. It is conve- 
nient (see Fig. 1) to represent a field (a tp field) by a leg with (without) 
an arrow. The interaction term ^int is a sum of three-leg vertices of types 
gi and §2 in Fig. 1. The nth order cumulant in Eq. ( p.5| ) comes from dia- 
grams with n such vertices. An external leg of a diagram will stand for a 
< field and an internal line for the contraction of two > fields. Each vertex 
is integrated over space and time. Only one-particle irreducible diagrams 
need to be considered. Moreover, any diagram in which one can go around 
a loop in the direction of the arrows vanishes by causality. A contraction of 
any pair of fields at equal times is zero, a property conveniently expressed 
by the convention 0(0) = and commented upon e.g. in Ref. Q. Finally, 
products of noncontracted fields at different times are reduced to single-time 
expressions by Taylor expansion in the time difference. 

After the integration on the momentum shell fields is done one applies to 
a scale transformation which depends on exponents z, d^, and d;^ still 
to be determined, 

x = bx' ip<{x,t) =b-'^^'ip'{x',t') 

t = bH' ^^{x,t) = b-'^^^'{x',t') 

and drops the primes of the new variables. The result is a renormalized 
action whose coupling constants are given in terms of the original couplings. 
In dimension d = dc — £ the cumulant series in Eq. ( |2.5D usually corresponds 
to an expansion in powers of e and upon truncating it at the desired order 
one has a RG involving only a finite number of coupling constants. If a fixed 
point of the transformation can be found, then scale invariance is a direct 
consequence and scaling laws follow. 



2.2 Stationary state of Siyp 

The time reversal symmetry of Eq. ( |1.18 ) imposes that d^ = d:^. We write 
d^ = = I (d + 77) , where rj is the anomalous dimension of the fields ip and 
ip. We have applied the renormalization procedure of subsection 2T to 5dp, 
taking into account all diagrams with two and three vertices. The resulting 
recursion relations read, with 6 = e^, 



dg 

di 



2 + 



Y 



2g' 



,^4 
A2" 



(2.7) 
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for the coupling constant of the three-leg vertex, and 

- = (2.8) 

for the deviation from the critical density. Here is the surface area of 
the unit sphere in divided by (27r)^. Requiring that the coefficients of 
the terms with the time and space derivatives in 5*0 remain constant leads 
to the conditions 

-^-5^:^ = (2.9) 



2A2 

^-2-r/-<7'^ = (2.10) 



These allow one to express the exponents z and r\ in terms of g^A 2. Eqs. (2.7- 
2.10| ) have a nontrivial fixed point, 

where we neglect terms of higher order in e. This fixed point is stable since 
in its vicinity g scales with the negative exponent yg = — e. The fixed point 
value of a is 

a* = ^ (2.12) 

In view of the definition of a this implies that the critical density pc in 
dimension 4 — e is higher than the mean field value p™^ = l/kr. Near the 
fixed point a and time scale with the exponents 

y^ = „-^ = 2-^ (2.13) 

z = 2-^ (2.14) 

respectively, and the anomalous dimension of the fields ijj and tp is 

r/ = -| (2.15) 

In Eq. (|1|) u is the usual correlation length exponent. Since /3 = ^^{d+r]), 
the density of the B species scales with 

13 = 1-^ (2.16) 





The analysis has been extended to second order in e by Janssen [|11| . 

The logarithmic correction in dimension d = 4, 

Pb{oo) aln^ a (2-17) 
has not, to our knowledge, been mentioned before. 
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2.3 Relaxation to the stationary state of Directed Percola- 
tion 



Relaxation to the stationary state of the directed percolation problem has 
not been considered before but can be studied by the same techniques. One 
has to keep the initial time term, with " coupling constant" , in the action 



of Eq. ( 1.2C| ). This term represents an initial state with Poissonian distri- 



bution of the B particles. In the diagrams it will be shown as an arrowed 
leg starting from a t = vertex indicated by an empty circle (see Fig. 1). 
It is easy to extend the analysis to a non-Poissonian initial distribution of 



particles. In that case we have to include in the action S dp of Eq. ( |1.20| ) 
another initial time term, viz. 

J d'^x V^^(x,0) (2.18) 

Higher order terms in the field ^/J also appear, but are irrelevant under renor- 
malization. The coupling constant in Eq. ( 2.1 j ) is related to the width 
of the distribution by 

= -^[(An^) - im)] (2.19) 

which vanishes for a Poisson distribution. The naive scaling dimension of 
the A^ term is easily found to be zero. This means that it is marginal and 
has to be included in the RG calculation. Any other terms generated by 
these two initial terms are irrelevant. Fig. 2 shows the two diagrams that 
yield renormalization contributions to /?^^ and to A^ which are linear in 
these quantities. The conclusions of this section may be arrived at on the 
basis of only these diagrams. Nevertheless we shall complete the calculation 
to one-loop order by also computing the contributions of the diagrams of 
Fig. 3. The diagram of Fig. 3(a) is the most complicated one to consider and 
the only one that we shall discuss explicitly. Its expression is 



A2 J{2^YJ ' ' ' VQ^(27r) ^2.20) 

One expands the fields in the integral for small t and r, and the exponential 
in powers of k. One can then carry the time integrations out explicitly and 
use the initial time constraint |Q] 

^{k, f = 0) = pg^ - 2A^V'(fe, t = 0) (2.21) 
It appears that out of all terms only 

2A^^pf J d^x ^{x, 0) (2.22) 



renormalizes pj^ and 



4A2^ / d<^x/(x,0) (2.23) 
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renormalizes A^. The other terms are irrelevant. The recursion relations 
for the couplings and that one finally finds are, to one loop order, 

^Pb^ 1/, ^ (o)^„. K^A^ (0) g'^K4 

- ^M-ri)py +2^^g -2p;^^A^— (2.24) 



2^ '^^^ ' '^^2A(A2 + ,t) "'^ A2 

^ = A,(-. - 2,^§ - VA,§ - A,,^§) (2.25) 



The first term in Eq. ( |2.24| ) and the term proportional to r] in Eq. ( p.25| ) 
are due to rescaling of the fields and the variables. The remaining terms in 
Eq. ( |2.24| ) come from the diagrams of Figs. 2a and 3a, and those in Eq. ( |2.25| ) 
from the diagrams of Figs. 2b, 3a, and 3b, respectively. Upon combining 
Eq. (|]2|) with Eqs. (|2Tll ) and (|1|) one sees that only the fixed point 



A^ = is stable, and A^ scales with the negative exponent j/A ^ = ~ f and 
hence is irrelevant. Substitution of this fixed point value in Eq. ( ^.24| ) shows 
that p'j^ scales with the exponent 

yB = \{d-7^) (2.26) 

— 2 

We note in passing that, since the same diagrams renormalize and 

^2, 



A^'i/' |^_Q, comparison of Eqs. (2.25) and ( |2.7| ) shows that the identity 



yA, = 2-z-| + | (2.27) 
must hold to every order in e. 

2.4 Initial and long time scaling of the relaxation 

We analyze the relaxation towards the stationary state following Diehl and 
Ritschel |14] (see also [^]). To simplify the discussion we place ourselves at 



the critical density, that is, put a = 0. Since psi't) is equal to the average 
of ip^Xjt) with respect to e~^°^, it scales with the exponent {d + rj)/2 and 
one has the scaling behavior 

PB{t) = b~^Fib-% b^p^^^) (2.28) 

valid in the limit 6 ^ oo and for fixed arguments of the function F. Choosing 
b^ = t and employing the relation /? = ^{d + r])/2 leads to 

p^(t)=t-^F(l,pg)til+^') (2.29) 

In the limit of large time i, at fixed initial density one obtains the usual 



critical relaxation proportional to t^^/'^^ at the condition that F{x, y) ~ 1 for 
y — > oo. In the limit of small at fixed time t a scaling regime appears 
whose existence was first pointed out by Janssen et al. 0]. In this limit 
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the density pB{t) should be proportional to p^^\ which is possible only if 
F{x, y) ~ y in the small y limit. Hence we can write 

PB{t)=P^^YHp^^h^^^') (2.30) 
in which J-{y) = y~^F{l, y) and 



e' = -"1 = — 
z 12' 

where the function T has the limit behavior 



(2.31) 



^(y)~i (y^o), Hy)-y'^ (y^oo) (2.32) 

Here 0' is the critical initial slip exponent Q. The fact that it is positive 
here implies that in an initial time regime (excluding microscopically small 
times) the B particle density increases as 

PB{t) ~ pff (2.33) 

Crossover to the regime of asymptotically long times takes place at t ~ Tcross, 
with 

W ~ ipf)-"^^^^''^ (2.34) 

The positivity of 9' and the initial rise of psit) have been discussed and 
interpreted in Refs. They can be qualitatively understood by the 

lack of critical fluctuations in the initial state. In the case of the reaction- 
diffusion systems considered here the fluctuations in the stationary state are 
responsible for the shift of the critical density to a larger value compared its 
mean field value. This is reflected by the negative cr* obtained in Eq. ( 2.12]) . 



Since the density distribution in the initial state is Poissonian and thus free 
of long range correlations the system shows qualitatively a mean field-like 
behaviour during the initial stage of the relaxation. For a = cj* < the 
density of B-particles therefore increases at (macroscopically) short times t. 
Notice, however, that this qualitative argument applies generally not to 



microscopically short times [16]. 



An important point, stressed in Ref. |Q] and not apparent from the pre- 
ceding discussion, is that 9' is an independent critical exponent. Whereas 
in the special case of ^dp considered here it could be expressed as —rj/z 
(see Eq. ( p. 31] )), that relation does not hold in general. This will be seen, in 
particular, in Section |3[ 

The scaling function T is known at the mean-field level, 

T{y) = ^ (2.35) 

1 + constant x y 

In dimension d = 4 one finds that the B density decays as 

pB{t) ~ Pb^ InS t {t < Tcross) 

PB{t) ~ t-Hnh (t> Tcross) (2.36) 
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where Tcross ~ (Pb^) Mna 

Finally, when the system is slightly off criticality, that is, for a ^ 0, the 
scaling function depends on the supplementary scaling variable |(T|t^/'^^. 
The power law decay ~ t~^/'^^ is then cut off exponentially on a time scale 
of order where ^ ~ (t~^ is the correlation length. 

3 Wilson renormalization of the full action S 
3.1 Introduction 

The RG procedure applied above to the action ^dp of the Directed Perco- 
lation problem can be extended without difficulty to the action S given by 
Eq. ( 1.16| ) and describing the full problem defined in Section 1.1. This action 



contains a parameter n which distinguishes between the diffusion constants 
of the healthy and the sick individuals. It appears to constitute an imme- 
diate and natural step beyond the action of Directed Percolation: Led by 
an interpretation different from ours, and starting from a pair of Langevin 
equations postulated for a problem in population dynamics, Kree, Schaub 
and Schmittmann Q arrive at the same action but specialized to /i = 0. Us- 
ing field-theoretic techniques they analyze its stationary state and compute 
the critical exponents 7y, and v to order e. 

Below we first derive the stationary state RG equations for the general 
case /i 7^ 0. Since p appears to be relevant, a new analysis is needed. For 
^ < we find a new fixed point and calculate the new critical exponents to 
lowest order in e. 

Next we derive the RG equations for the additional couplings that de- 
scribe the initial state. We shall, in particular, consider an initial particle 
distribution which is independent from site to site, but non-Poissonian. 



3.2 Stationary state of the full action S 

Our starting point is the action S of Eq. ( |1.16D . Since our first interest 
is in the stationary state, we ignore for the moment the initial time term 
proportional to . As before we write S = Sq + S'lnt- The free part of 



the action now includes the term in Eq. (1.16) that is proportional to p, and 



whose presence is due to the diffusion constants of the A and B particles 



being unequal. This term violates the time reversal symmetry of Eq. ( 1.18| ) 
and introduces a nonzero correlator 

(A - IjA;^ + Xa 

Diagrammatically this contraction is represented by a solid leg with an arrow 
that connects to a dashed leg without arrow to form an internal line (see 
Fig. 1). As a consequence the symmetry between the V'^V' ^ind V'V' vertices 
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in the action is broken and we write the interaction part Si^t as 



'int 



-I 



,^2 



(3.2) 



with independent gi and 52 • 

Some further general considerations are useful. The vertices of S\-at (see 
Fig. 1) do not allow for renormalization of the (p propagator. We therefore 
conclude that 



z = 2 



(3.3) 



In order to show that this is the dynamic scaling exponent not only of the 
conserved density ip but also of the fluctuations of the order parameter, one 
has to check that the coupling A tends to finite nonzero fixed point value 
(see below). From dimensional analysis the dimensions of ip and Tp must be 
equal to |. However, because of the symmetry breaking discussed above, 
one has to allow for two distinct anomalous dimensions, r] and fy, of the fields 
■0 and 

We have worked out the RG recursion relations for the couplings, taking 
into account all one loop diagrams (they involve up to three vertices) . Some 
of these diagrams are shown in Fig. 4. For convenience we put 



X 



9192 



K4, Y 



92UH 



4A2 (A + 1)2 ^' A2(A + 1)2^ 

With this notation we are led to the following recursion relations. 



d£ 



--91 



£_^_|_8X + (3 + i)y-(5A + 3)Z 
(A + 1)3 ^2 ^ {\ + lfYZ 



X 



4A X 



(3.4) 



(3.5) 



dff2 
d£ 



92 



£_|_^_8X + (3 + ^)y-2(2A+l)z' 



(3.6) 



dn 
dZ 



= u 



\-\{v + f])-^X + Y-{2\ + l)z] 



(3.7) 



1, , 1 A2-4A-1 



(3.8) 



Upon requiring that the the coefficient of the f/'^ti/' term remain equal to 
unity we find the supplementary equation 



= -i(^ + fy) - 2X + y - i(3A + l)Z 



(3.9) 
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Finally, the coefficient Act renormalizes according to 



dl ' 2^' 2A A2 + (j 

- ^ TTT-^^rrTTZ (3-10) 



+ 



(A + 1)A2 + Act 
A (A2 + a)[(A + l)A2 + Ac7] 



The preceding equations, together with the appropriate initial conditions for 
£ = 0, determine the RG trajectories. By combining Eqs. ( ^.10 ) and 



we may check, before even solving the full set of equations, that to order e 

= 2 - I - 4X* + - (2A + = 2 - I (3.11) 

is independent of fi. The existence of a Ward identity (see the appendix in 
[^) which is preserved for fJ- ^ guarantees in fact that this result is true 
to all orders in e. 

3.2.1 Equal diffusion constants, ^ = 

For ^ = the time reversal symmetry is restored, we must set r] = fj, and 
the RG equations Eqs. (|3.5D and ( |3.6D for gi and g2 become identical. The 
critical behavior is controlled by the fixed point 

X* = |, y^ = |, Z^ = 0, A = 2, r/ = -| (3.12) 

first found by Kree et al. Q. 

3.2.2 Unequal diffusion constants, fJ, 



The flow defined by Eqs. ( |3.5| - |3.8D and subject to condition (|3.9|) possesses 
a stable fixed point at which A is the solution of 

2A^ - OA^ - 6A - 1 = (3.13) 

that is 

A = [(2 + Vsy/^ + (2 - ^3)^/^ - 2]-^ = 5.10664... (3.14) 
and furthermore 

X* = '4±^e = 0.227e (3.15) 

8A(3A + 1) ^ ' 



Y* = = 0.418e (3.16) 
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A — 2A — 1 
^'- 2A(A+l)(3A + l) -' = °'°"°-' 

This leads to the critical exponents 

71 = 77 = — e= -0.313 e (3.18) 

' ' 3A + 1 ^ ^ 

The result rj = holds at every order in e. It follows that the exponent 
P = v{d + ri)/2 takes the mean field value (3 = 1. In checking the linear 
stability of this fixed point in the subspace of {gi, g2,u, X) we found two 
negative eigenvalues of multiplicity two. It is easy to see that this fixed 
point governs the critical behavior of only the stationary state with fi < 0. 
Indeed, the unrenormalized {i = 0) value of Z has the sign of — and from 



the definition of Z combined with Eqs. (3^) and (|3.7| ) one deduces that Z 
cannot change sign under renormalization; so that the fixed point can be 
reached only if Z{£ = 0) is positive just like Z*. 



3.3 A first-order transition for > ? 

We shall now assume that the stationary state of the system undergoes a 
first order transition for /i > 0. This assumption makes of /i = a tricritical 
point. We now develop a scaling analysis of the stationary state for fj, of 
either sign. In this subsection the exponents describing the ^ < fixed point 
will bear a minus index. 

At the tricritical point ^ = characterized by r], u, and (3 the coupling 
constant /i is a relevant perturbation with crossover exponent = —1^/2 > 
0. For both a and p small compared to microscopic momentum scales pB 
therefore has the scaling form 

PB{a,p) = {-afn^^{-a)-'^y-) (3.19) 

The critical behavior for o" — > at fixed /.i < is described by the exponents 
rj- and /3_ and we get 

Pb{oo) ~ A{-af- (3.20) 

The amplitude A is nonuniversal and is a singular function of p, because 
the limits cr — > and — > do not commute. If p is small compared to 
microscopic scales but large compared to {—ay^f^ we see from Eq. ( p.igj ) 
that 

A ~ (-;^)(/5-/5-)/(^3/m) (3.21) 

We use analogous arguments in the case /i > assuming that the phase 
transition is first order. A first order transition means that pb{oo) has a 
finite value in the limit —a — > at fixed ^ > 0. Eq. ( |3.19| ) now gives 

Pb{oo) ~ (3.22) 

with 5 = yy^/ [3 = —i]/{d + rj) > 0. Fig. 6 shows a qualitative plot of the 
resulting behavior of the stationary state density pb{oo) as a function of p 
in the three regimes p < 0, p = 0, and ^ > 0. 
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3.4 Relaxation to the stationary state of the full action S 

3.4.1 Non-Poissonian initial states: general remarks 

We wish to consider an initial particle distribution which, as before, is the 
product on i of identical single-site distributions p{mi,ni), where rm and rii 
count the A^s and i?'s, respectively, on site i. The distribution p is arbitrary. 
It need not be Poissonian in either rrii or n^, and need not factorize. The 
interest of studying such initial distributions is that all initial distributions 
with short-range correlations reduce to them under renormalization. The 
general distribution of this class is represented in the action S by the initial 
time terms 

d'^x[-pgV(a:, 0) + A^Tp'^ix, 0) + A^^ifix, 0)'^{x, 0) + A^V^^(x, 0)] 

(3.23) 

plus terms of higher order in Tp and ■0, which are irrelevant under renormal- 
ization. Here, setting ii = nii + rii, we have 

= -^[(An?)-(n,)] 

= -^KA£f) 

A<^^ = --[{An^i) - {n,) + {AmiAm)] (3.24) 
P 

If the A and B particles have uncorrelated Poisson distributions, then A^p = 
A^pip = A^ = 0. 

3.4.2 Non-Poissonian initial A distribution 

For the sake of simplicity we focus on the fi = case. As discussed 
above, in an arbitrary initial state we have to deal with the four param- 
eters p^^\ Aip, A^^, and A^. Also for simplicity, and because it is one of 
the interesting cases, we shall consider here A^ = A^^ = 0. Eq. ( |3.24| ) 
implies that this corresponds to having a Poisson distribution of the B par- 
ticles, a non-Poissonian distribution of the A particles, and no correlations 
between the occupation numbers of the two species. It then follows that 
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A^ = - — [(Am2)-(m,)] (3.25) 



The A^ term allows for contractions of the p between themselves according 
to 

{p{h,h)p{k2,t2)) = -A^(27r)'^5('^)(fci + A:2)e(ti)e(t2)e-'^?(*i+*^) 

(3.26) 

There are no diagrams by which the A^ term gets itself normalized: it is 
strictly marginal. However, the new contractions of Eq. ( p. 26 ) contribute to 
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the renormalization of the term. The additional diagrams that appear 
are shown in Fig. 5. The exphcit expression of diagram Fig. 5(b) reads 

U'A^ f / dtdri^ik, t)lp{-k, t + T)e{t)e{T) f ^-2<l't-X[(,^kf+a]r~^r 

One repeats the analysis performed in subsection |2.3| ; however, unlike dia- 
gram Fig. 4(a) that appears in the directed percolation case, diagram Fig. 5(b) 

yields a renormalization of that does not vanish at the fixed point. The 

recursion relation for the vertex is 



where we have set 



rjo = {X + l)Y''A^ = -A^e (3.29) 

o 

From Eqs. and (|]2|) we conclude that the scaling dimension of p^j^ 

is 

This exponent is therefore nonuniversal. We may repeat the scaling argu- 
ment of section 2A, which allows us to conclude that the critical initial slip 
exponent is 

n, d + rj 7] r/o 
& = ^ Vr = 

= 3.31 

32 ^ ' 

where in the last equality we have used that z = 2. We have therefore 
exhibited a nonuniversal exponent that characterizes the short-time behavior 
of the density of the B species. From expression ( 3.25| ), combined with the 
fact that at criticality the density of A particles differs negligibly from p, 
one may deduce that 

1 

2' 



< -P (3.32) 



so that 9' has the (nonuniversal) lower bound ^^j^ = (4 — 9/o)e/64. The 
coupling A<^ attains the maximum value allowed by Eq. ( ^3.321) , and hence 9' 
its minimum, when the rrii have the probability law 

PAim) = {I- p)Sm,,0 + pSm„l (3.33) 

(where we take p < 1). This means that the A particles are randomly 
distributed but without any multiple occupancy. Apparently, this is the 
environment giving the slowest initial rise of the number of infected individ- 
uals. 
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3.4.3 Relaxation for < 



We shall now consider critical relaxation to the stationary state for /i < 0, 
but limit ourselves to Poissonian initial distributions: = = A^^ = 0. 
The initial B particle density renormalizes according to 



A (0) 1 

^ = l(d-v)P^' (3-34) 



After slightly generalizing the reasoning of subsection 2.4 we conclude that 
there is a nontrivial critical initial slip exponent 

To summarize, in the initial stage for ^ < 0, psit) ~ P~b^^'^ then, over a 
time scale the density decays as pB{t) ~ t^'^^'^. This last result holds at 
every order in e. 



4 Conclusion 

We have adapted Wilson's renormalization scheme to the two-species reaction- 
diffusion process A + B ^ 2B, B ^ A. The stationary state of this process 
exhibits a second order phase transition at which the S's become extinct. 
A parameter /i = 1 — Db/Da controls the relative strength of the diffu- 
sion constants and Db of the two species. We have determined the 
critical behavior of the stationary state for < 0, that is, when the the 
species subject to extinction (B) diffuses faster than the species (A) which 
it parasites on. Special emphasis has been put on the relaxation from var- 
ious initial states to the critical stationary state. As happens in magnetic 
systems, we found an initial relaxation regime which scales with an inde- 
pendent exponent 9' , arising from renormalization of the initial time terms 
in the action. The reaction-diffusion process has an order parameter (the 
B particle density) that couples to a conserved quantity (the total density). 
We have shown that 0' varies continuously with the width of the distribution 
of the conserved density. The same effect occurs in magnetic systems with a 
continuous symmetry of the order parameter, where the symmetry implies 
the existence of a conserved quantity coupling to the order parameter . 

We found that for positive ^ the renormalization transformation has 
no stable fixed point, which signals the absence of a second order phase 
transtion. We therefore discussed the possibility of a first order transition; 
this scenario is plausible but we lack decisive arguments in its favor. This 
unresolved point deserves further attention. 

A question of great interest concerns the low-dimensional behavior of 
this reaction-diffusion system. The e expansion is meaningful slightly be- 
low the upper critical dimension dc = 4, and one may imagine that it still 
gives a fairly good idea of what happens in d = 3. However, the qualita- 
tive picture provided by the e expansion may break down below some new 
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critical dimension d'^. Eg. ( 1-19| ) shows that the (unrenormalized) action 
contains terms, neglected in this work, that become relevant with respect 
to the free (Gaussian) theory when d = 2. It is therefore possible that low 
dimensions require separate consideration |l^. Even the seemingly simple 
one-dimensional version of this model exhibits nontrivial behavior |19| that 
remains to be elucidated. 
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FIGURE CAPTIONS 



Figure 1. Graphical conventions. The tp propagator is shown in plain line, 
the if propagator appears as a dashed line, and whenever it is non-zero, the 
{ip(p) correlator connects a plain leg to a dashed one. Also shown are the 
vertices comprising the interaction term of the full action. A vertex con- 
necting legs by an empty circle o is attached to t = and is therefore not 
summed over time. 

Figure 2. Diagrams (a) and (b) renormalize and A^, respectively. 
They are linear in these quantities. 

Figure 3. Other one- loop order diagrams renormalizing p^'* (diagram (a)) 
and (diagrams (a) and (b)). 

Figure 4. Additional diagrams entering the renormalization of gi for 7^ 
with respect to the p = case. 

Figure 5. Diagrams (a) and (c) renormalize Ay,. Diagram (b) renormalizes 
and A^. 

Figure 6. Stationary state density pb(oo) of the B particles plotted quali- 
tatively against the total density p. The curve for /it > (dashed line) rests 
on the hypothesis of there being a first order transition. The curve for fi = 
(dotted line) was obtained in Ref.[l]. The curve for /x < (plain line) is a 
mean-field like straight line. 
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